Analysis of bone erosions in rheumatoid arthritis using HR-pQCT: Development of a measurement algorithm and assessment of longitudinal changes

Purpose The purpose of this study was to establish an algorithm for measuring bone erosions at metacarpophalangeal (MCP) joints using high-resolution peripheral quantitative computed tomography (HR-pQCT), to investigate the precision of measurements, and to assess longitudinal changes in bone erosions among patients with rheumatoid arthritis (RA). Methods The 2nd and 3rd MCP joints were scanned at a voxel size of 60.7 μm using second-generation HR-pQCT. Bone erosions on MCP joints were identified using a semi-automated algorithm we developed, and each erosion parameter was measured. Measurement reproducibility was evaluated in 19 healthy subjects using intraclass correlation coefficients (ICCs) and root mean square percent coefficient of variance (RMS%CV). Finally, longitudinal changes in bone erosions over a period of 12 months were assessed in 26 patients with RA based on the calculated least significant change (LSC). Results Reproducibilities for measurement parameters regarding bone erosions with our algorithm were good (all ICCs ≥ 0.98; all RMS%CVs < 5%). No erosion parameters showed significant changes after 12 months of treatment in terms of median values in all erosions, while both progression and repair of erosions were observed individually (e.g., erosion volume: progression, 26% (+0.62 mm3); repair, 34% (-0.85 mm3); no change, 40%). Conclusions The measurement algorithm developed for bone erosions at MCP joints showed good reproducibility. Both progression and repair of bone erosions were observed in patients with RA even after 12 months of appropriate treatment. Our algorithm may be useful to investigate the etiology of RA and assess drug efficacy.


Introduction
Rheumatoid arthritis (RA) is a systemic chronic inflammatory disease that causes cartilage and bone destruction in association with polyarthritis and leads to disability [1]. Bone erosion, as one of the characteristics of RA, is an important predictor for severity of disease and prognosis related to bone destruction [2,3]. The identification of bone erosions in early RA and precise assessment on subsequent follow-ups are essential for diagnosis and treatment.
In recent clinical practice, computed radiography (CR) has been widely used to assess bone erosions [4]. Magnetic resonance imaging (MRI) and ultrasonography (US) are more effective for identifying bone erosions than CR [5]. However, CR only provides two-dimensional (2D) information, and MRI and US have limitations in providing detailed three-dimensional (3D) evaluations of bone.
High-resolution peripheral quantitative computed tomography (HR-pQCT) allows in vivo analysis of bone microarchitecture at peripheral sites with high spatial resolution [6,7]. Currently, this modality is used for various studies in bone diseases, particularly in the field of osteoporosis [8,9].
Research regarding RA using HR-pQCT has been increasing since around 2010 [10,11]. In particular, HR-pQCT is advantageous for identifying small bone erosions compared to CR and MRI, because of the higher resolution [12]. Although bone erosions were clearly identified using HR-pQCT for not only RA patients, but also healthy subjects, indicating that bone erosions are not specific to RA [11,13], bone erosions remain an important indicator in RA patients. Studies using HR-pQCT focusing on the analysis of longitudinal changes (progression and repair) in bone erosions are increasingly being conducted to verify the effects of medical treatments for RA [14][15][16][17][18][19][20][21]. On the other hand, several analytical methods have already been reported for identifying bone erosions and evaluating measurement parameters [22][23][24][25]. However, the precisions of these methods have not been fully established.
The main purposes of this study were to: 1) establish methods for analyzing bone erosions in metacarpophalangeal (MCP) joints using second-generation HR-pQCT; 2) verify the precision of parameters for measuring bone erosions using this algorithm; and 3) assess longitudinal changes in bone erosions over 12 months in RA patients applying this algorithm.

Subjects
Forty-eight healthy women with no history of finger arthritis were recruited to verify the precision of our algorithm in this study. Another 26 patients with RA fulfilling the 2010 American College of Rheumatology/European League Against Rheumatism (ACR/EULAR) classification criteria [26] participated in the longitudinal study and received treatments continuously with the same medications determined by a rheumatologist according to EULAR recommendations during the observation period.

Scans
The 2nd and 3rd MCP joints of one hand were scanned using a second-generation HR-pQCT system (XtremeCT II; Scanco Medical, Bruttisellen, Switzerland) with reference to the scan protocol recommended by The Study group for xtrEme Computed Tomography in RA (SPECTRA) [27]. In healthy subjects, the same hand was scanned three times with repositioning on the same day, to verify the precision of our algorithm. In patients with RA, the affected hand was scanned at baseline and at the 12-month follow-up. Image quality was evaluated according to the motion artifact grading [28], and images with grade 4-5 were excluded from subsequent analysis.

Identification of bone erosions
In this study, erosions at the MCH and phalangeal base (PB) of the 2nd and 3rd MCP joints were analyzed using bone microstructure measurement software (TRI/3D-BON; Ratoc System Engineering, Tokyo, Japan).
An overview of our algorithm for identifying and analyzing bone erosions is shown in Fig  1. First, an image of the MCH or PB to be analyzed was selected ( Fig 1A) and a binarized image was created using a threshold value of 320 mg/cm 3 , which is a standard value for trabecular bone in second-generation HR-pQCT ( Fig 1B). The bone mask image was then automatically created by filling the bone marrow space (Fig 1C). Subsequently, all concavities on the bone surface were automatically filled (Fig 1D), and the contour of this image was shrunk by 0.25 mm (Fig 1E). All concavities on the bone surface were then extracted by subtracting Fig  1C from Fig 1E (Fig 1F). As shown in Fig 2, the following concavities were excluded as falsepositive patterns: physiological concavities (Fig 2A); pseudo-concavities formed by osteophytes ( Fig 2B); and small cortical interruptions or vascular channels ( Fig 2C). Finally, the remaining concavity was determined to represent a bone erosion after checking on a 2D axial view (Fig 1G), a coronal or sagittal view (Fig 1H), and a 3D view (Fig 1J and 1K) with reference to the definitions of erosion defined by SPECTRA [22].

Assessment of identified bone erosions
The volume of bone erosion was automatically calculated in voxel units based on the extracted volumes of interest (VOIs) (Fig 1L). Widths and depths of bone erosions were manually measured based on 2D images of the VOI. The depth and width of bone erosion in the axial view were defined as Er.D and Er.W, respectively ( Fig 1G). In addition, depending on the distribution of erosions based on four quadrants, as published earlier [11], the depth and length of bone erosion in the perpendicular (coronal or sagittal) view were defined as Er.D-2 and Er.L, respectively ( Fig 1H). (i.e., bone erosion on the radial or ulnar side was measured in the coronal view, and bone erosion on the dorsal or palmar side was measured in the sagittal view). In addition, peripheral bone mineral density surrounding bone erosion (Er.BMD) was automatically measured based on the VOI extracted with a width of 1 mm (Fig 1I) [24].

Additional procedures
Additional procedures are necessary for some erosions, as shown in Fig 3. First, the step shown in Fig 3A may be necessary if overflow of the VOI beyond the concavity occurs (Fig 3Aa and 3Ab). In this case, the appropriate VOI was extracted by shrinking the contour of the image corresponding to Fig 1E by 0.50 mm (Fig 3Ac and 3Ad), which was applied because values above 0.50 mm made the VOI too small.
Second, the steps shown in Fig 3B may be necessary if the VOI connects with bone marrow spaces (Fig 3Ba). In this case, the VOI needs to be temporarily and repeatedly reduced in size until the continuities are completely interrupted (Fig 3Bb) [24]. Dilation operations are then executed until the VOI returns to its original shape (Fig 3Bc). Finally, the newly selected VOI   is defined as the VOI of a bone erosion (Fig 3Bd). To ensure measurement reproducibility, the same number of size-reduction operations in these processes is used in longitudinal analyses.
Third, the adjustment of rotation in the axial view may be necessary if a certain degree of difference in rotation exists in the longitudinal analysis. As shown in Fig 3C, even using the same images, different rotation causes slight differences in filling concavities (Fig 3Ca and  3Cc) and the extracted VOIs (Fig 3Cb and 3Cd).

Statistical analysis
Intra-reader reliability for each bone erosion parameter using our algorithm was evaluated using the intraclass correlation coefficient (ICC) and root mean square percent coefficient of variance (RMS%CV) [29]. Least significant change (LSC) was calculated by multiplying the intra-reader precision error by 2.77, and both percent LSC (LSC of percent coefficient of variance: LSC%CV) and absolute LSC (LSC of standard deviation: LSCSD) were evaluated [30].
The Wilcoxson signed-rank test was used to assess longitudinal changes of bone erosions between baseline and 12 months follow-up in patients with RA. Er.BMD was evaluated only for bone erosions identified at baseline.
In addition, cumulative probability plots were created to analyze longitudinal changes of individual bone erosions in patients with RA. LSCSD calculated for each parameter was applied to these plots, and the percentages of the three categories of "progression", "repair" and "no change" were calculated based on each LSCSD we regard as cut-off values. Similarly, percentages of these three categories (increasing, decreasing, and no change) for Er.BMD were calculated based on LSCSD for this parameter.
All statistical analyses were performed using JMP version 14 (SAS Institute Inc., Cary, NC, USA) and values of p<0.05 were considered significant.

Measurement reproducibility of bone erosions
Twenty-five erosions were identified in 19 of 48 healthy subjects using our algorithm, and 15 erosions showing motion artifact of grade 1-3 in three scans in a row (Fig 4) were used to verify the reproducibility of measurement parameters. These results are shown in Table 1. ICCs of all parameters showed almost perfect reproducibility (0.98-1.00). RMS%CV of each parameter was less than 5% (1.91-4.65%), indicating good reproducibility. LSC%CV was 5.3% for Er. BMD, and 9.98-12.87% for all other parameters. LSCSD was 0.21 mm 3 for Er.V and 14.7 mg/ cm 3 for Er.BMD, and 0.13-0.18 mm for width and depth of bone erosions.

Characteristics of patients with RA
Characteristics of patients with RA (20 women, 6 men) at baseline are shown in Table 2. Mean age was 66.5 ± 9.2 years, and disease duration was 108.9 ± 143.6 months. The numbers of patients divided into disease durations of <12 months, 12-120 months, and >120 months were 8, 10, and 8, respectively. Mean height was 1.56 ± 0.08 m, weight was 54.5 ± 10.9 kg, and BMI was 22.3 ± 3.9 kg/m 2 . RF-positive results were seen in 88% of patients, and ACPA-positive results in 50%. Mean parameters at baseline and follow-up were 0.69 and 0.11 mg/L for CRP, 33.9 and 15.9 mm/h for ESR, and 3.7 and 1.7 for DAS28-CRP, respectively. Methotrexate was being administered to 54% of patients, and biologics to 81%. In addition, 27% of patients were using low-dose glucocorticoids.

Longitudinal changes of bone erosions in patients with RA
Longitudinal changes of bone erosions between baseline and follow-up (12 months) in patients with RA are shown in Table 3. Thirty-four bone erosions (MCH: 23; PB: 11) were confirmed at baseline, and 37 (MCH: 26; PB: 11) were confirmed at follow-up, including 4 new erosions and complete repair of 1 erosion. For the total cohort, no significant changes between baseline and follow-up were seen for any parameter.
Longitudinal changes of each bone erosion are shown in Table 3

Discussion
In this study, MCP joints were scanned using second-generation HR-pQCT (voxel size, 61 μm), and a semi-automated analysis algorithm was established for bone erosions. In addition, the precision of each parameter regarding bone erosions was verified, revealing good Table 1. Measurement reproducibility of each erosion parameter.

Mean (SD) Median (IQR) ICCs RMS%CV (%) LSC%CV (%) LSCSD
Er.V (mm 3  Three scans and measurements were repeated for the reproducibility analysis. Values are calculated for precision error and displayed as mean and standard deviation reproducibility (all ICCs �0.98; all RMS%CVs <5%). Finally, longitudinal changes of bone erosions after 12 months in patients with RA were investigated. Each parameter showed no significant changes in the total cohort, while both progression and repair of erosions were observed on an individual despite appropriate treatment. Because HR-pQCT enables detailed analysis of bone erosions, changes of bone erosion volume (progression and repair) are the most commonly used indicators for disease activity or treatment efficacy in RA. However, methods for analyzing bone erosions currently vary between different research groups.
One method was to estimate erosion volume from an ellipsoid formula based on manual measurements of depth and width for a visually identified erosion. This approach has been in wide use from the beginning of RA research using HR-pQCT [22,23]. In addition, semi-automated algorithms to improve measurement precision have been reported, such as the 3D segmentation technique described by Töpfer et al. [24], or the cortical interruption detection described by Peters et al. [25]. In terms of intra-reader reliability for assessing bone erosion volume, RMS%CV was 5.66% for the algorithm described by Töpfer et al. [24], and ICC was 1.00 for the algorithm described by Peters et al., although RMS%CV was not mentioned [25]. In terms of the measurement reproducibility of bone erosions in the present study, RMS% CV was 3.97% and ICC was 1.00 in terms of the intra-reader reliability of our algorithm for assessing bone erosion volume. These results were comparable with those for the other two algorithms. One of the strengths for our algorithm is the measurement of depth and width based on the VOI. Although seemingly the same as the manual measurement mentioned above, the VOI clearly delineates sites of bone erosion in our algorithm. The RMS%CV of depth and width was thus less than 5%, indicating high precision.
Regarding longitudinal changes of bone erosions in this study, no parameter for bone erosions showed significant changes between baseline and follow-up for the total cohort, while both progression and repair of erosions were observed at the individual level. In particular, the results in this study showed progression in 26%, repair in 34%, and no change in 40% for erosion volume.
Existing  Values at baseline and follow-up, and differences between them are displayed as median and interquartile range (IQR).
The Wilcoxon signed-rank test was used for comparisons of each value between baseline and follow-up (significant p < 0.05). progression, and no change in volume of interruptions were 7%, 6% and 87%, respectively [20]. Yue et al. reported that a reduction in erosion volume was apparent after 1 year, and the proportions of patients with erosions showing progression, partial repair and no change were 32%, 20% and 48%, respectively [32]. Such findings may have resulted from differences in study designs, such as background characteristics including duration of disease or pharmacotherapies, along with different methods of analysis. Nevertheless, a certain proportion of longitudinal changes in volume showing both progression and repair in each erosion was observed even with appropriate treatment, and our study supported the findings of those previous studies. We expect that the algorithm described here will contribute to pathological analyses of RA or the effects of specific pharmacotherapies. Some studies have already verified the effects of specific pharmaceutical agents such as denosumab and tocilizumab [17,21], and our algorithm may contribute to similar research in the future.
This study has several limitations. First, we verified the precision of bone erosion measurement in healthy subjects, not in RA patients, due to the need to scan several regions such as the fingers, wrist, radius, and tibia on the same day, in consideration of the problems of scan time and radiation dose. Bone erosions in healthy subjects were considered similar to those in patients with RA, as shown in Fig 4, and have been used as controls in previous studies [11,20]. Second, we were unable to analyze huge bone erosions connected to bone marrow spaces. We were also unable to identify and analyze erosions that were too small. Small erosions not identified by our algorithm were considered as minute vascular channels or cortical interruptions, not bone erosions. Finally, images with severe motion artifacts may require careful analysis because of the potential for insufficient contouring of the periosteal line, particularly in longitudinal studies.
In conclusion, we established an analytical algorithm for bone erosions on MCP joints using second-generation HR-pQCT. Good reproducibility was obtained with our algorithm for measurement parameters regarding bone erosions. We observed both progression and repair of bone erosions at the individual level with longitudinal follow-up after 12 months in patients with RA, despite appropriate treatment. We believe that our algorithm allows precise assessment of bone erosions and may be useful for investigating the pathogenesis of RA and the effects of treatment.